One-dimensional transport of interacting particles: Currents, density profiles, phase 

diagrams and symmetries 
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Driven lattice gases serve as canonical models for investigating collective transport phenomena 
and properties of non-equilibrium steady states (NESS). Here we study one-dimensional transport 
with nearest-neighbor interactions both in closed bulk systems and in open channels coupled to two 
particle reservoirs at the ends of the channel. For the widely employed Glauber rates we derive an 
exact current-density relation in the bulk for unidirectional hopping. An approach based on time- 
dependent density functional theory provides a good description of the kinetics. For open systems, 
the system-reservoir couplings are shown to have a striking influence on boundary-induced phase 
diagrams. The role of particle-hole symmetry is discussed and its consequence on the topology of 
the phase diagrams. It is furthermore demonstrated that systems with weak bias can be mapped 
onto systems with unidirectional hopping. 

PACS numbers: 05.50.-|-q, 05.60.Cd, 05.70.Ln 



I. INTRODUCTION 

One-dimensional driven transport has manifold appli- 
cations in biology, physics, and materials science. Promi- 
nent examples are the motion of motor proteins along 
microtubules or actin tracks [l|, y] , protein synthesis by 
ribosomes [3|, ion diffusion in narrow channels [^-la] or 
charge transfer in photovoltaic devices [3, Q- Many of 
those have been studied by models based on incoherent 
hopping processes, where either the focus was on an ef- 
fective one-particle description Q, or on the collective 
behavior of mutually excluding particles as describe d by 
the asymmetric simple exclusion process (ASEP) pi]| - [l3l |. 

From the fundamental point of view, one-dimensional 
driven systems are of vital interest also to gain a better 
understanding of the physics of non-equilibrium steady 
states (NESS). These are macrostates carrying steady 
currents. An important question is whether and how 
concepts and theorems well known for equilibrium sys- 
tems, for example, the fluctuation-dissipation theorem, 
Onsager reciprocal relations, and maximum entropy con- 
siderations can be generalized to NESS (l4| - [lq . Most 
challenging is certainly the question, if, as in equilibrium 
systems, a limited number of control variables can be in- 
troduced, which allows one to make general statements 
with respect to the distribution of microstatcs or struc- 
tural and kinetic properties of NESS. As kind of "minimal 
models" , totally asymmetric simple exclusion processes 
(TASEPs), where particles can hop only in one direction, 
are particularly suited for corresponding studies. 

The standard TASEP refers to particles on a one- 
dimensional lattice, which mutually exclude each other 
and perform jumps to vacant nearest neighbor sites to 
the right with a rate F. Considering a bulk system, 
easily realized by employing periodic boundary condi- 
tions, where A^p particles occupy on a large ring of N 
sites, corresponding to a density p = Np/N. In this case 
one finds that the distribution of microstates is uniform, 
that means all particle configurations are equally prob- 



able [13, [l3 ■ The bulk current j^ is thus exactly given 
by the mean-field expression jg = F(A^p/A^)[l — (A^p — 
1)/{N - 1)], yielding Jb = Tp{l - p) in the hmit of infi- 
nite system size. The process becomes more interesting 
when considering an open system, where particles are in- 
jected from a left reservoir with particle density pi^ and 
ejected to a right reservoir with particle density p^. In 
this situation the distribution of microstates in the NESS 
is no longer uniform but can be calculated analytically 
by utilizing a special matrix algebra lid. 131. recursion 
relations [lE, ll9|, or the Bethe ansatz |lC)l.ll2|. 



Moreover, there as an intriguing phenomenon appear- 
ing in the NESS, namely the bulk density pe far from 
the boundaries to the reservoirs shows phase transitions 
as a function of the control variables Pl and pu_. Phase 
diagrams can be derived from so-called minimum and 
maximum current principles |20l . |2l|. These principles 
state that if the particle density Pl is lower (higher) than 
the density p^, a bulk density ps is established in the 
system, which corresponds to the minimum (maximum) 
of the bulk current Jb(p) in the range Pl < p < Pr 
(Pr < p < Pl). They are a consequence of the fact that to 
match the reservoir densities at the boundaries, density 
profiles in the system cannot be uniform in general, and 
accordingly changes in the bulk current must be com- 
pensated by diffusive currents. For example, if Pl > Pr, 
and the local density is assumed to decrease monotoni- 
cally from the left to the right, then the diffusive current 
should be positive everywhere, or zero in regions of con- 
stant density. Accordingly, the current in the bulk region 
of flat density proflle must be at a local maximum. It is 
important to realize that this argument relies on the as- 
sumption that the density profile varies monotonically. 

In further studies j2ll - |23j it has been shown that the 
minimum and maximum current principles can also be 
applied to certain TASEPs with particle-particle interac- 
tions beyond (athermal) site exclusions. In some analogy 
to equilibrium systems, this suggests that the bulk behav- 
ior is determined by experimentally controllable reser- 



voir properties and independent of microscopic details of 
system-reservoir couplings. However, as was shown re- 
cently |24l |. application of the minimum and maximum 
current principles requires a very specific way of particle 
injection and ejection in this case. 

In general, density oscillations appear at the system 
boundaries in the presence of interparticle interactions, 
which implies that the minimum and maximum current 
principles can no longer be used to predict boundary- 
induced phase diagrams |25| . To capture the transport 
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in the presence of such oscillations requires a theory that 
allows one to connect correlations to the density profile 
on a local scale. The time-dependent density functional 
theory (TDFT) of lattice gases [H, [H is well suited for 
this situation. In particular, combined with the Markov 
chain approach to derive microstate distributions in equi- 
librium as functionals of the density [2^, [29[ it allows 
one, in a rather straightforward manner, to calculate re- 
lations between correlators and densities in equilibrium 
systems with inhomogeneous density profiles. As a con- 
sequence, the method becomes a powerful means to de- 
scribe kinetics, and we will refer to it as the Markov chain 
approach to kinetics (MCAK) in the following. It has 
the merit that it becomes exact for bulk kinetics with a 
canonical distribution of microstates. Using the MCAK, 
boundary-induced phase diagrams can be predicted with 
good accuracy. The resulting phase diagrams appear to 
be very different for different system-reservoir couplings, 
not only with respect to locations of transition lines but 
also with respect to the overall topology. This finding 
is somewhat surprising, in particular because it seems at 
first glance that particle-hole symmetry gets broken for 
nearest-neighbor interactions. One of the goals of this 
work is to clarify the reason for the change in topology 
and the associated question regarding particle-hole sym- 
metry. 

A further goal is to study whether the results reported 
in [231 for TASEPs with interactions remain valid for 
ASEPs, where jumps against the bias direction are pos- 
sible, as it is the case in any realistic application. In this 
connection we also reanalyze the driven transport when it 
is mediated by Glauber jump rates, which, among other, 
have been used in the field of incoherent electron trans- 
port along molecular wires [3y, |3l|. Interestingly, for 
these Glauber rates an exact expression can be derived 
for the bulk current-density relation. This is because 
the Glauber rates belong to a class, where a canonical 
Boltzmann distribution is valid for the microstates in the 
NESS. We also demonstrate that the MCAK not only 
provides good descriptions of the NESS but also of the 
dynamic time evolution of density profiles. 



II. TASEP WITH NEAREST-NEIGHBOR 
INTERACTIONS 

We consider a one-dimensional lattice gas with hard- 
core exclusion, unidirectional nearest-neighbor hopping 
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FIG. 1. (Color online) Sketch of particle jump with rate 
r(ni_i, ni+2) from site i to site {i + 1) for the different possi- 
bilities of occupations of sites (i — 1) and (i + 2). 



with rates Fi^j+i and repulsive nearest-neighbor interac- 
tion V > 0. The microstate of the system is specified by 
the set of occupation numbers n = {rii}, where each site 
i of the system is either occupied by a particle (n^ = 1) 
or vacant (n^ = 0). The total energy of the system is 
given by the lattice gas Hamiltonian 



n = vj: 



niUi+i . 



(1) 



Using the master equation for the time evolution of the 
probability density P{n,t) of microstates, the evolution 
equations for mean values Pi(t) = (ni)t = X^n "'^^('^'^) 
(henceforth called densities) are j32| 



dt 



ii-l.i{t) - Jhl+lit) : 



(2) 



where ji^i+i{t) is the average current from i to {i + 1), 
jt,i+i{t) = (?ij(l - ni+i)Ti.i+i{n))t . (3) 



Here, 



refers to an average over P{n,t). As il- 



lustrated in Fig. [H the rates Fi_i_|_i(n) are functions of 
the occupation numbers ni_i and ni+2 only, Ti^i+i{n) = 
F(ni_i, 711+2)- Accordingly, the current in ([3]) can be 
written explicitly in terms of four-point correlators, 

jt,i+i = (?lj-i7ijri.i+inj+2>tr(0,0) 
+ (7ij_ini7ij+irij+2>tr(l,0) 
+ {ni-inih,+in,+2)t r(0, 1) 
-f (nj_ininj+inj+2)tr(l,l) . (4) 

Here we introduced hole occupation numbers fii ~ 1 — n^. 
As mentioned in the Introduction, we here use the 
widely employed Glauber rates [33j 
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where v is an attempt frequency, /3 is the inverse thermal 
energy, and AH = {ni+2^ni^i)V is the energy difference 
between the states after and before the jump. In the 
following we set /3 = 1 and v ~ 1. Because F(0, 0) = 
F(l, 1), the bulk dynamics is particle-hole symmetric, i.e. 
a bulk system with particle concentration p and the set of 
jump rates {F(0,0), F(1,0),F(0, 1), F(l, 1)} is equivalent 
to a bulk system with particle concentration 1 — p and 
the set of jump rates {F(l, 1), F(l, 0), F(0, 1), F(0, 0)}. 



III. BULK CURRENT-DENSITY RELATION 

To evaluate the bulk current-density relation in NESS, 
one has to determine the correlators in Eq. Q . In general 
this is a difficult task because, different from equilibrium 
systems, there are no universal laws yielding the distri- 
butions of microstates in NESS. On the other hand, some 
authors [3J, [S^ have considered the question whether it 
is possible to specify the rates r(ni_i, 71.^+2) in such a way 
that the distribution of microstates in the NESS equals 
the equilibrium Boltzmann distribution ex exp(— "H). In- 
deed it was found that this is the case, if the rates satisfy 
the relations 



r(o,i) = r(i,o)e-^, 

r(o,o) + r(i,i)-r(o,i)-r(i,o) 



0. 



(6a) 
(6b) 



A derivation of these relations is given in the Appendix, 
because it was not given in detail in the original work 

Interestingly, the Glauber rates satisfy Eqs. (pal) and 
(j6bp . As a consequence, the correlators in Eq. Q equal 
the equilibrium correlators in the corresponding one- 
dimensional Ising model, which can be calculated by vari- 
ous means, such as the transfer matrix technique, density 
functional theory etc. The result for the current reads 
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(7) 



where / = l/[exp{V)+l] = r(0, 1) and C^^'> = {n,n,+i),^ 
is the equilibrium nearest-neighbor correlator. 
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(8) 



Because the bulk dynamics is particle-hole symmetric, 

J(p)=j(l-P)- 

Figure [5] shows the behavior of the current as a function 
of density for various interaction strength V . For y ^> 0, 
j{p) approaches the parabola j = {p — P^)/2 for parti- 
cles with site exclusion only. When V exceeds a critical 
value 14 = 21n3 ~ 2.20, j(p) develops a double-hump 
structure 1201 l36i with two maxima at densities 



PhiV) = 2 T 



\ 




(9) 



and a minimum at half-filling, i.e., for p = 1/2. In the 
limit case V — > 00, we find j = {x^^-^ — 2x + x^^^)/{2 — 2x) 
with X = (2p — 1)^, meaning that there is no particle 
movement for p = 1/2. For pl^ = 1/2 T (\/2 - l) /2 the 
current is maximal, in agreement with earlier findings 
reported by Krug |20| . 



0.14 
0.12 
0.1 
0.08 
0.06 
0.04 
0.02 






^^-^^^ 








^'—--^ 


\ 




/ 




v = 


o\ 


- / 


\ • — 


V = 


^* \ 




\ / — ■■ 


V = 


2K \ 




■.: 


V^oo \ 








\ 



0.2 



0.4 0.6 

P 



0.8 



FIG. 2. (Color online) Bulk current-density relation j{p) for 
various interaction strengths V. 



IV. TRANSPORT IN OPEN SYSTEM: 
APPLICATION OF MCAK 

Coupling of the TASEP to a left and right reservoir 
in general requires eight coupling parameters for nearest- 
neighbor interactions, as indicated in Fig.|31 For injection 
of particles to site i ~ I, uq and ai specify the injection 
rates if site i = 2 is vacant or occupied, respectively. Due 
to the missing neighbor on the left for particles on site i = 
1, in addition the rates for the two possible jumps from 
site i = 1 need to be specified. These are denoted by a^ i 
for vacant/occupied site i = 3. Analogously, /3o,i denote 
the two possible ejection rates for vacant/occupied site 
i = N — I, and /3q ^ the two possible rates from site 
i = iV — 1 for vacant/occupied site i = N — 2. 

To evaluate the currents in Eq. ^ , wc cannot use any 
longer the mapping of the NESS to an (unbiased) equi- 
librium state as discussed in the previous Sec. IIIIl be- 
cause the translational invariance used in the derivation 
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FIG. 3. Couplings of the system to the (a) left and (b) right 
reservoirs mediated by the a and j3 rates. Indices and 1 
refer to the occupation of the sites next to the target site (for 
a rates) and to the initial site (for j3 rates) . The primed rates 
are for jumps from and to boundary sites of the system. 



(cf. Appendix) is broken. As known also from the stan- 
dard TASEP with site exclusion only, the distribution 
of microstates is not uniform in the open systems, that 
means it changes when going from the bulk to the open 
system. 

To treat the relevant correlators in Eq. (|4]) one can 
consider their time evolutions. This would lead to the 
appearance of higher-order correlators and different pro- 
cedures could be applied for closing the resulting hi- 
erarchy. However, this approach usually becomes un- 
handy Instead we use the underlying concept of TDFT 
[2a,[23,l33l , which is based on the (time-)local equilibrium 
approximation. This amounts to approximate the non- 
equilibrium distribution P{n,t) by the Boltzmann prob- 
ability ex exp[—'H{n)] plus an effective time-dependent 
external potential J^i hi[p{t)]ni, where p{t) = {pi}. This 
implies that the correlators at any time t are supposed 
to be related to densities as in an equilibrium system. 
These relations are now needed for inhomogeneous sys- 
tems without translational invariance. In particular, as 
mentioned in the Introduction, it is important to include 
information on the local variation of the density. 

To this end the Markov approach for express- 
ing the distribution of microstates [28| is particularly 
suited. In this approach the equilibrium joint prob- 
abilities peq [rii, ■ . ■ ,ni+j) for the occupation num- 
bers rii, . . . ,ni^j are expressed by the Markov chain 

piq{ni,...,ni+j) = Poq{ni)Yll^^w{ni+s\ni+s^i), 
where poq {rii) is the probability for rii in equilibrium, and 
w{nij^i\ni) = Pcqini,n^+i)/p[,f{{ni) is the conditional 
probability for n^+i given rii. Since the joint probabilities 
are directly connected to the correlators, e.g., piq (fii-i = 
0, n,; = l,ni+i = 0,ni+2 = 0) = {ni^inini+ini+2)cq, aU 
four-point correlators in Eq. ^ can thus be reduced to 
two-point correlators. 

Applying this MCAK procedure to the TASEP with 
nearest-neighbor interactions yields 



Selecting the physical branch of the solution, we obtain 
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+ ct\c^l\rio,i) + cl'_\cll\T{i,i] 

where 

C> = {nini+i)cq ^ Pi- C> , 
Q = {ninij,i)cq = Pi+i - Cf' , 

Q = (^i"-i-|-l)oq ^l - Pi- Pi+1 + Cf' , 

and C follows from the quadratic equation 
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(11a) 
(lib) 
(lie) 

(12) 



c 



(1) _ 



1 



2(1 -e-y) 



(p,+p,+i)(l-e-^)-l 



+ \/[(p. + P»+i)(l - e-y) - if + Ap,p,+ie~y{l - e~y) 

(13) 

In addition to the currents ji,i+i not directly coupled 
to the injection and ejection rate, we need the currents 
at the boundary sites. 



.ks == {nih2)aa + {fiin2)ai , 
iia = ("-ifi2"3)ao + ('^i^2"-3)ai , 
JN.iK. = ("-Af-i"Af)/3o + {nN-inN)Pl , 
JN-1,N = ("-Af-2«Af-l"-Af)/3o 

+ {nN-2nN-inN)l3'i ■ 
Using the method outlined above, we obtain 

Jl,i = CI 'ao + Cl 'ai, 

/^(2) 

3i,2 = C'2 'a^ + C'2 '«! , 

i — /92 ^ / 

JAT.R = Cj^_^I3q 4- Cj^_^Pi , 



C 



JN-1,N 



(2) 



PN- 



~ \pN-2Pii + ^N-2l^lj 



(14a) 
(14b) 
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Given the explicit expressions PU)) . P^ for the currents 
in terms of the densities via Eqs. (|TT|) . ((T2|) the kinetic 
equations ^ become a closed set. 



V. BOUNDARY-INDUCED NESS PHASES 

As mentioned in the Introduction, applicability of the 
minimum and maximum current principles requires spe- 
cific "bulk-adapted couplings" of the system to the reser- 
voirs. In fact, the a and /3 rates need to be defined in such 
a way that the system can be viewed as being continued 
into the reservoirs, corresponding to relations between 
correlators and densities as in the bulk. 

In a bulk system, when an initial configuration 
{ui^i = 0, ni+2} would be given, two rates are 
possible for a particle jump from site i (i.e. rii = 
1): r,,j+i = l/[exp(n,+2F) + 1], if n,_i = 0, 
while Fi^j+i = l/[exp(r7,j+2^ - V) + 1], if nj_i = 
1. For given {rii+i = 0,71.^+2}, let us denote 
by p(01|0nj+2;p) = p(010ni+2;yo)/p(0ni+2; p) and 
p{ll\Qni+2; p) = p(110ni+2;p)/p(0ni+2;/3) the condi- 
tional probabilities for the configurations {ni_i,ni} = 
{0,1} and {ui^i^rii} = {1,1} to occur in the NESS 
of a closed bulk system with density p and interaction 
V. For example, injection rate ao.i then results from 
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FIG. 4. (Color online) Boundary-induced phase diagrams of NESS at V = 2V* for (a) the bulk-adapted and (b) the equilibrated- 
bath couplings. Thick solid and thick dashed lines mark first- and second-order phase transitions for the TASEP. Corresponding 
thin lines are for the ASEP with F = 2. In (b) the symbols refer to KMC results and the lines to MCAK results. 



a weighting of rates with the probabilities p{01\0n2', pi^) 
and p(ll|0n2; Pl) corresponding to virtual configurations 
{ti-i = 0, TiQ = 1, rii = 0, 712} and {?i_i = 1, uq = 1, 711 = 
0,712} at the boundaries. Following the same procedure 
for the other rates we arrive at (777 = or 1) 

p(01|077i;pO , p(11|0777;pO 



exp{mV) + 1 exp[(77i — 1)V] + 1 

p{0\lOm;p^) ^ p{l\10m;p^) 
exp{mV) + 1 exp[(77i — 1)V] + 1 



p(00|l777;p:,) 



p(10|l777;pn) 



exp(-777y) + 1 exp[(l - m)V] + 1 



p(0|0l777;p,,) 



p(1|01777;pr) 



exp(-777y) + 1 exp[(l - m)V] + 1 



(16a) 



(16b) 



(16c) 



(16d) 



Here, p(0|107n;p) and p{l\10m; p) are, respectively, the 
bulk probabilities for 77i_i = and 77i_i — 1 under the 
condition that {77^, 77i_|_i, 77^+2} = {1,0, 7(7,}. Note that 
for the P rates the given occupation numbers are those 
to the left side, i.e. p(00|l777;p) = p(77il00;p)/p(7r7l;/9) 
and so on. 

Application of the minimum and maximum current 
principles to the TASEP with V = 2\4 and the bulk- 
adapted rates in Eqs. (J16p yields the boundary-induced 
phase diagram shown in Fig.UJa). In total seven phases 
occur, where the bulk density equals either the left reser- 
voir density Pl (phases I and VI), or the right reservoir 
density p„ (phases III and V) , or the densities p^ 2 [see 
Eq. (|n])] of maxima in the current (phases II and VII), 
or the density 0.5 of the (local) minimum in the current 
(phase IV). Transitions between these phases can be of 
first or second order, which are indicated by thick solid 
and thick dashed lines, respectively. The thin lines re- 
fer to changes of the phase diagram when allowing for 



jumps against the bias direction, as further discussed in 
Sec. IVIII Notice that the diagram has symmetry with 
respect to the diagonal p^ = 1 — p^, which reflects the 
particle-hole symmetry in the system as explained in the 
following Sec. I VII 

The results for the bulk densities and currents of the 
NESS in the case of bulk-adapted couplings and for our 
choice of rates satisfying the relations ([6]) are exact. Note 
that this does not hold true for the density profile close 
to the boundaries. Application of the MCAK described 
in Sec. IIVI allows one also to calculate the time evolution 
of density profiles. Corresponding numerical solutions 
of Eq. ^ with the expressions for the currents derived 
in Sec. IIVI are approximate both at the boundaries and 
in the bulk, because at transient times the relations be- 
tween correlators and densities differ from those in the 
equilibrium state without bias. 

In order to get insight how well the MCAK captures 
the kinetics, we have performed kinetic Monte Carlo 
(KMC) simulations of the TASEP with bulk-adapted 
couplings for a chain of A^ = 1000 sites with Pl = 0.9, 
Pr = 0.6, V = 2T4 and an initially empty lattice. Re- 
sults from these KMC simulations (symbols) for density 
profiles at five different times, as well as the stationary 
state, are compared in Fig. [SJa) with the predictions of 
the MCAK (lines). Three different time regimes can 
be distinguished. The first regime is the "penetration 
regime" for t < 2000 during which the initially injected 
particles pass the system and reach the right reservoir. In 
this regime there is excellent agreement of the KMC data 
with the MCAK predictions. The penetration regime is 
followed by an "intermediate regime" , where the density 
in the system increases until approaching values close to 
the hmiting one in the NESS. The two times t = 4000 
and t = 8000 in Fig. [SJa) belong to this regime. For the 
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FIG. 5. (Color online) Time evolution of pi for (a) bulk-adapted and (b) equilibrated-bath couplings of the systems to the 
reservoirs for /9l = 0.9, pn = 0.6, V = 2Vi,, N = 1000, and pi{t = 0) = 0. KMC results are marked by symbols and MCAK 
results by lines. In the KMC simulations, averages were performed over 10^ different configurations. The inset in (b) zooms 
out the oscillatory behavior at the boundaries. 



choice of parameters in the present example, it is interest- 
ing that a kind of domain wall appears in the system, see 
the jump-like change of density at i « 900 for t = 4000 
that has moved to i « 200 for t = 8000. A second such 
kind of domain wall appears in the time interval 8000 
and 16000 and moves to the right (not shown). 

One can view the occurrence of these transient domain 
walls as resembling the occurrence of domain walls along 
first-order lines in the phase diagrams of the NESS. Con- 
trary to the latter, the positions of the transient walls not 
only fluctuate, but they exhibit an average drift, because 
the local current in the system is not constant. For the 
wall seen in Fig. [SJa), the current left to the wall must 
on average be larger than right to the wall. The MCAK 
captures the formation of transient domain walls, but 
the quantitative agreement with the KMC data is less 
accurate than in the penetration regime. The interme- 
diate regime is followed by a "relaxation regime" , where 
at each point the density continuously relaxes, without 
rapid jump- like changes, towards the limiting value in the 
NESS. In this regime the MCAK predictions are again in 
excellent agreement with the KMC data. 

The bulk-adapted couplings are specifically tuned to 
make the minimum and maximum current principles ap- 
plicable. With respect to applications such couplings 
will not be realized, but one is led by the fact that 
the time scale of relaxation processes in the reservoirs 
is much faster than in the system. With this assump- 
tion, baths can be assumed to correspond to equilibrated 
Fermi gases with chemical potentials /Xl = ln[pL/(l — Pl)] 
and /Zr = ln[pf,/(l — Pn)]- A reasonable ansatz for the a 
and (3 rates then is 



a„ 



/9L[exp(my - /11) + 1] ^ 
[exp{mV) + 1]"S 



(17a) 
(17b) 



/3™ = (1 - PR)[exp(/iR - mV) + l]-^ , (17c) 

P'^ = [eM~mV) + 1]-^ (17d) 

The Fermi factors in these rates correspond to the 
Glauber rates, if one considers that injected particles 
loose an energy ^l, ejected particles gain an energy /Ir, 
and that the interaction with particles in the system is 
as in the bulk. The additional factors Pl in Eq. (|17ap 
and 1 — Pr in in Eq. (|17cp take into account the filling 
of the baths. The functional form in Eq. ([T7| resembles 
forms resulting from Fermi's golden rule for transition 
rates [30, ISlI, lla 123 ■ We will refer to the couplings me- 
diated by the rates in Eq. pT)) as the "equilibrated-bath 
couplings" . 

For these couplings, the density profiles at the bound- 
aries can no longer be expected to vary monotonically, as 
it is required for applicability of the minimum and max- 
imum current principles. Considering equilibrium sys- 
tems, it is well known that modified interactions, for ex- 
ample at confining walls, commonly lead to density oscil- 
lations. It would be surprising if such density oscillations 
do not appear for NESS under modified interactions at 
the boundaries, as for the equilibrated-bath couplings. 

Figure [5]Jb) shows the time evolution of density pro- 
files obtained from KMC simulations (symbols) and the 
MCAK (lines) for the same reservoir densities and cou- 
pling strength as in Fig.ISJa). Indeed, density oscillations 
appear at the walls in the stationary state, as demon- 
strated in the inset of Fig. JUb). They can be under- 
stood when pointing out that the oscillations are miss- 
ing in the bulk due to translational invariance. If one 
would, in the bulk, determine the spatial dependence of 
the density with respect to an occupied site, which in 
fact amounts to a determination of density correlations, 
then it is clear that oscillations occur due to the repulsive 
nearest-neighbor interactions. The same holds true when 



the spatial dependence of density profiles is determined 
by starting from a vacant site. The reservoir in case of 
the equilibrated-bath couplings resembles a vacant site 
(missing nearest neighbors) and therefore oscillations ap- 
pear at the boundary. In agreement with this picture, the 
density at site N next to the right reservoir is particu- 
larly large, the density at the next site iV — 1 to the left 
then particularly low, and this alternating behavior con- 
tinues on the scale of the correlation length towards the 
bulk. The bulk density is p^ — 0.34 and deviates from 
the bulk density Pb — 0.70 in Fig.ISJa) following from the 
maximum current principle. 

To cope with the oscillations a theory is needed where 
the local current ji,i+i is dependent on the form of the 
density profile around sites i and i -I- 1 , as it is the case in 
the MCAK. As can be seen from Fig. ^h), the density 
oscillations at the boundaries are well accounted for by 
the MCAK, and accordingly the predicted bulk density 
is in excellent agreement with that from the KMC sim- 
ulations. Also the time evolution of the density profiles 
is well captured by the MCAK in Fig. [SJb). Let us note 
that the breakdown of the minimum and maximum cur- 
rent principles does not imply that phases corresponding 
to the local minimum and to the maxima in the bulk 
current density relation can no longer appear. In fact, 
starting from the flat region and considering the onset of 
the bending of the profile at the ends of this region, an 
enlarged region of monotonically varying density profile 
could be considered, where the minimum and maximum 
current principles apply. 

The boundary-induced phase diagram for V = 2V^ and 
the equilibrated bath couplings is displayed in Fig.|4l[b). 
This diagram strongly differs from the corresponding one 
for the bulk-adapted couplings in Fig. HJ^a) . Instead of 
seven phases, five phases appear, where the bulk density 
either is determined by the left reservoir density Pl via 
a function /i(pl) (phase I), or is determined by the right 
reservoir density pR, via functions /iii(pr) and /v(Pr) 
(phases III and V), or is equal to p* (phase II), or is 
equal to the density 0.5 of the (local) minimum in the 
current (phase IV). Analogous to Fig. IHa), first- and 
second-order transitions are marked by thick solid and 
thick dashed lines, respectively. Phases in Fig. ^a) and 
(b) labelled by the same Roman numbers correspond to 
each other in the sense that their character agrees (left or 
boundary determined, or minimum or maximum current 
phases). In addition, the phase diagram in Fig. |4l[b) is 
not symmetric with respect to the diagonal Pr = 1 — 
Pl. The reason for this will be clarified in the following 
Sec.lVD 



VI. PARTICLE-HOLE SYMMETRY 

Under exchange of particles by holes the ejection rates 
/3o, /3i, /3o, P'l would correspond to injection rates ai, 
ao , ct'i J ct'o I respectively, and the current direction would 
be reversed. Because the bulk dynamics is particle-hole 



symmetric, the particle-hole exchanged system must have 
the same properties with respect to the hole occupation 
numbers fii = 1 — ti^, i.e. densities Pi{t) = {ni)t and 
correlators {nifijhk ■ ■ ■)t at any time t in the particle-hole 
exchanged system equal pi{t) = {ni)t and {nkn-jTii . . .)t in 
the original system. In this sense particle-hole symmetry 
holds true in general. 

The bulk density in particular must fulfill in the NESS 

Psiao, ai, u'q, a[, l3o, l3i, P'q, l3[) (18) 

= 1 - Pb (/3i , /3o , /3i , /3o , ai , ao , a'l , "o) • 

For compact notation, let us introduce multivariate a 
and 13 jump rates under particle-hole exchange of the 
a = (ao,ai, aQ,a[) and /? = (ao, ai, aQ,a[) rates, 

«=(ai,«o,«'i,ao), /3=(/3i,/3o,/?;,/3^). (19) 

Then we can rewrite Eq. ([T8)) as pB(a,/3) = p-aiP^a) = 
1 -pB(/3,a). 

Following the view that the reservoirs are controlled by 
only a few variables, as their chemical potentials or den- 
sities, we should require the injection and ejection rates 
to depend on Pl and Pr, respectively. Given a = ^(pl) 
and /3 = /3(pr), the bulk density becomes a function of 
Pl and Pr, 



/5b(Pl,Pr) = pB(a(pL),/3(pR)) 



(20) 



Particle-hole symmetry would show up in this function, 
if the relation 



/5b(Pl,Pr) = 1 - Pb(1 - Pr, 1 - Pl 



(21) 



is fulfilled. Replacing the left hand side with Pb(Pl, Pr) = 
PB(a(pL),^(pR)) = l-pB(/3(pR),a(pL)) by using Eq. HH]), 
and the right hand side by 1 — Pb(1 — Pr, 1 — Pl) = 1 ^ 
PB(a(l — Pr),/3(1 — Pl)), we obtain by comparison 



aip) = f3ip) 



(22) 



as the condition for the particle-hole symmetry in 
Eq. (ini) to be obeyed. 

The rates in Eqs. (|17p for the equilibrated-bath 
couplings do not satisfy relation (P^ and hence the 
phase diagram in Fig. |H[b) does not display the 
symmetry according to Eq. (HI]). For the bulk- 
adapted couplings, in contrast, the rates in Eq. (J16p 
satisfy Eq. (|^ . For example, /3i(l — p) = 

(exp(-\^)+l)-ip(1100; l-p)/p(ll; l-p)+2-^p{110l; 1- 
p)_/p{n;l- p) = {cM-V) + l)-ip(1100;p)/p(00;p) + 
2 ^p(0100; p)/p(00; p) = ao, where we have used the 
particle-hole symmetry p(ni_i,n,i,rii+i, nj+2; 1 — p) = 
p{fii-i,fii,hi+i,fiij^2', p) in the bulk. Analogously, the 
other relations in Eq. ([^^ can be proven. 

Theoretically, given the particle-hole symmetric bulk 
dynamics, the behavior in the open system is fully con- 
trolled by the eight a and /3 jump rates, and the different 
boundary-induced phases would appear in a particle-hole 



symmetric m.anner in the respective eight-dimensional 
space. However, in practice it will be difhcult to "control" 
couplings in this detailed way. Rather the system will be 
connected somehow to the reservoirs and one could tune 
the reservoir properties. In the case considered here, this 
is reflected by Eq. ([20| . which parameterizes the eight 
rates in terms of two densities. As a consequence, differ- 
ent phases from the eight-dimensional space are projected 
out into the (pl , yOR)-plane for different coupling mecha- 
nisms. This can go along with significant changes of the 
topology, as indeed obtained in Fig. S) 



VII. ASEPS 

TASEPs are simplified models because jumps against 
the bias direction are not included. For hard-core exclu- 
sions only, it is known that the structure of the boundary- 
induced phase remains essentially the same when allow- 
ing for backward jumps [40, |41| . In the presence of a bias 
F in forward direction, forward and backward jump rates 
r_i. and r.j_ are generally assumed to fulfill the detailed 
balance condition, i.e. F^/F^ = cxp(— AH), where the 
lattice gas Hamiltonian in the presence of the bias reads 



as well as the superscripts (2) and (3) in the correlators. 



H = V^^ i-iiTii+i - f y^ii 



(23) 



Considering a corresponding ASEP in the limit _F — > oo, 
an associated TASEP with rates F = limF_,.oo F^ is ob- 
tained, if the forward rates saturate for infinite bias. Con- 
versely, given a TASEP with rates F, an ASEP can be 
defined that in the limit F ^ oo reduces to the TASEP, 
for example, by setting F^ = F independent of F and 
F^ = Fexp(-A-H). 

Interestingly, TASEPs can even be associated with 
ASEPs in the linear response regime of weak bias F. 
Let us consider an ASEP with the Glauber rates from 
Eq. (O, where the forward jump rate is now given by 



r^(ni_i,ni+2) 



j(ni+2-ni-l)V-F _|_ Y 



(24) 



G 



(3) 



Ji+l,* = 



(cl:\cSF^(0,0) + Ca\caF^(0,l) 

+ cf^.CfXvui. 0) + cii\cl!\F^(i, 1) 



(25) 



Except that the ji 
rents 



_i must be replaced by the net cur- 



Jia+i{i) = ji,i+i{t) - ji+i,x(t) 



(26) 



between sites i and i + 1, the rate equations ^ remain 
the same. 

Bulk current-density relations in the NESS from the 
MCAK are compared to KMC results in Fig. |nia)-(c) for 
various interaction strengths V and three different bias 
values i^ = 1, 2 and 3. As expected, for small F = 1 
[Fig. [HJa)], the MCAK gives excellent agreement with 
the KMC simulations for all V. With increasing F, de- 
viations become significant for F >V, see, for example, 
the results for 1/ = K and F = 3 in Fig. Eljc). Note, 
however, that for the special case V = 0, correspond- 
ing to the standard ASEP with hard-core exclusion only, 
the MCAK always gives exact results, independent of F, 
because in this case all microstates are equally probable 
[1G |. In the regime of strong interactions V > F, the 
MCAK provides good results, and in particular agrees 
with the KMC results in the limit y — > oo. 

The current as a function of the bias F is shown in 
Fig. Eld) for one representative particle density pe — 0.3 
and three different values V = 0, Vi, and 214. For 
_F < 1, the current increases linearly with F, while for 
F > 1 nonlinear response effects become relevant. For 
large F, the currents from the KMC simulations satu- 
rate at values independent of V, while the limiting cur- 
rents in the MCAK are T^-dependent. The critical values 
Vc{F), where the bulk-current density relation develops a 
double- hump structure, sec Figs. [Sl^a)-(c), increase with 
stronger bias. In the inset of Fig. IHd), we display the 
MCAK results for Vc{F). Surprisingly, for F — >■ 0, where 
the MCAK becomes accurate, Vc{F) approaches the crit- 
ical value T4 for the TASEP considered in Sec. |lll 

To understand this, let us write the net current as 



and the backward jump rate by F^(ni_i, ni+2) = 
F_>(r7,i_i, ni+2) exp [(ni+2 — J^i-i)!^ — ^]- Because the 
F_i,(ni_i,ni+2) no longer satisfy Eqs. ([S]), calculations 
for the bulk behavior in the NESS, based on equilibrium 
relations between correlators and densities, are not exact. 
We can expect, however, that the MCAK will provide a 
good approximation for small bias F. 

The forward current ji.i+i from site i to site i -I- 1 has 
the same form as in Eq. ([T0| with the jump rates given 
by Eq. (|24p. The backward current ji+i^,; follows from in- 
terchanging Pi and Pi+i, and the indices (i — 1) and (j-l-l) 



J = p(0100)F^(0,0)- 


-p(ooio)r. 


-(0,0) 


-Fp(1100)F_ 


>(i,o)- 


-p(0011)F^ 


-(0,1) 


+ p{omi)T_ 


.(0,1)- 


-p(1010)F^ 


-(1,0) 


+ p(1101)F_ 


>(1,1)- 


-p(ioii)r^ 


-(1,1) 



(27) 



where we have combined in each line "reversed configu- 
rations" , i.e. equivalent situations for forward and back- 
ward jumps. The p(....) are independent of F in the 
MCAK (in the bulk) and the same for a given configura- 
tion and its reversed. Hence we can rewrite Eq. (|27p in a 
way that in each line the first p(....) are multiplied by the 
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FIG. 6. (Color online) (a)-(c) Current-density relations j(p) for the ASEP at various V and three different bias strengths F . 
In (d) the current-bias relation is shown for p — 0.3. Symbols refer to KMC data, lines to MCAK results, and the assignment 
of lines and symbols to the different interaction strength V is given in the legend of (a). 



difFerences r^(., .) — r<_(., .) between forward and back- 
ward rates. Taking the linear response limit of these dif- 
ferences gives J = Fjtasep + 0{F'^), where Jtasep refers 
to a TASEP with rates 



by Qf_j. and /3^ rates. The backward currents are 



r(ni_i,ni+2) = 



1 + cosh[(ni+2 - n^-i)V] 



(28) 



Again these rates do not satisfy Eq. (|6]), implying that 
the MCAK treatment of the bulk NESS behavior of 
this TASEP is no longer exact. Although the rates in 
Eqs. dMl) and (0 are different, the MCAK yields the 
same current-density relation given in Eq. ([7]). Accord- 
ingly, Vc{F) becomes T4 hi the limit F — >■ 0. 

That in the MCAK, the ASEP in the linear response 
regime can be associated with a TASEP raises the ques- 
tion whether this would be true in an exact treatment. 
Considering a general expansion of the right hand side of 
Eq. ((27|) for small F, this requires the p(....) to exhibit 
no linear terms in F. We have not yet achieved to prove 
this property, but representative KMC results shown in 
Fig. [7] for the p{--..) in the second line of Eq. (P7|) are in 
agreement with it. 

Let us now extend our discussion to open systems. The 
functional form of the boundary currents in bias direction 
are as in Eqs. (|15ap - (jl5dp with the a and (3 rates replaced 



Ji, 



J2,l 



Jn,N 



.(2) 



C['>a^,o + C{''a^,i, 



(1). 



CI 



(3) 



P2 



(Cf<0 + Cf^a^^,), 



Cjy_il3^,o + C^_^(ij, 



-(2) 



(^(3) 
jN.N-1 - I <^JV-2P^,0 + <-'Ar-2P^,l j • 



(29a) 
(29b) 
(29c) 
(29d) 



As discussed above, the MCAK provides an accurate 
description in the linear response regime, and the bulk 
behavior of the ASEP in this regime is equivalent to a 
TASEP with rates (gS]). Moreover, the MCAK predicts 
the same bulk-density relation for the rates psp as for 
the rates ([5]) . Accordingly, it is insightful to compare the 
boundary-induced phase diagrams of the ASEP with the 
TASEP in Sec. HH 

To this end we have, for the bulk-adapted coupling, ap- 
plied the minimum and maximum current principles to 
the bulk current-density relations of ASEPs for various 
F a.t V = 2V^, as, for example, to those shown Fig.[Sl^a)- 
(c). For the equilibrated-bath couplings we use Eqs. P7)) 
for the a^, /3^ rates with FL from Eq. (P51) and a<_. 
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FIG. 7. (Color online) KMC results for the joint probabilities 
p(llOO) andp(OOll) of finding a configuration {ni_i = l,ni = 
1, Wi+i = 0, ni+2 = 0} and its reversed configuration {ui-i — 
0, rii = 0, rii+i = l,ni+2 ~ 1} in the bulk part of the NESS 
for the ASEP defined in Eq. ((24)) . The interaction strength is 
V — 2T4 and the particle concentration p — 0.5. 



/3<_ rates determined by the detailed balance condition. 
The corresponding rate equations are integrated numer- 
ically for various F at V = 214 and the resulting density 
profiles analyzed in the long-time limit. 

For small F ^^ 0, we found that MCAK results for the 
phase diagrams in Figs.|3l^a) and (b) are almost the same 
for the ASEP. Visible small differences appear when F 
leaves the linear response regime. To illustrate this, we 
have indicated in both Figs. IH^a) and (b) phase transi- 
tions (thin solid and dashed lines) for F = 2. 



VIII. CONCLUSIONS 

Effects of interparticle interactions beyond hard-core 
exclusions in collective driven transport pose many chal- 
lenges and possibilities, whose significance has not yet 
fully explored. Here we have considered ASEPs and 
TASEPs with repulsive nearest-neighbor interactions in 
one dimension. For these jump processes on lattices, cur- 
rents can in general be expressed in terms of correlators of 
occupation numbers whose order increases with the inter- 
action range. To arrive at closed sets of kinetic equations, 
one has to decide on how to treat the relevant correla- 
tions. The Markov chain approach for deriving exact den- 
sity functionals [28| allows one to express correlators in 
terms of densities, where the respective relations, strictly 
valid in equilibrium, entail information on the local den- 
sity variation, which is necessary to capture interaction- 
induced non-monotonic behavior of density profiles. Let 
us note that a standard TDFT treatment based on an ex- 
act functional would provide such relations only via the 
solution of integral equations connecting the correlators 
with direct correlation functions. 

As we have demonstrated, application of the MCAK 
leads to a good description of both the time evolution 
of density profiles and their limiting shape in the NESS. 
Because of this, boundary-induced phase transitions of 



the bulk density as functions of reservoir densities could 
be well predicted. The coupling to the reservoirs turned 
out to have a decisive influence also on the topology of 
phase diagrams. Particle- hole symmetry in the nearest- 
neighbor interacting lattice gas with open boundaries 
manifests itself in certain relations between the injection 
and ejection rates. It was clarified under which condi- 
tions the particle-hole symmetry shows up also with re- 
spect to the reservoir densities in the boundary-induced 
phase diagrams. Furthermore we have demonstrated that 
ASEPs in the linear response regime can be mapped onto 
TASEPs with rates that are related to the first term in 
an expansion of the difference between forward and back- 
ward rates with respect to the bias. As a consequence, no 
significant changes in boundary-induced phase diagrams 
occur when connecting ASEPs with weak bias to corre- 
sponding TASEPs. 

For the jump rates we have used Glauber forms in this 
work. These were shown to belong to a class, where 
the distribution of microstatcs in the NESS is equal to 
the Boltzmann distribution of the interacting lattice gas 
without bias. Accordingly, an exact bulk current-density 
relation in NESS could be derived. We notice that the 
mapping of a NESS to an equilibrium state without 
bias would not be possible in higher dimensions for the 
Glauber rates |34| . 

In Refs. [13, [13 we considered TASEPs with jump 
rates ex exp(— AH/2), where, as in Eq. (O, AH is the 
energy difference between states after and before the 
jump. Different from the Glauber rates, these rates 
are not bounded and do not fulfill Eqs. ([Sa|) and (j6b[) . 
Nevertheless the phase diagrams for bulk-adapted and 
equilibrated-bath couplings are very similar to the ones 
displayed in Fig. [H This suggests that the bulk dynam- 
ics has only a weak influence in contrast to the dynamics 
coupled to the reservoirs. This suggestion is reinforced 
by the fact that the topology of the phase diagram ap- 
pears to be the same (for given boundary couplings), even 
if a bulk dynamics is considered that reflects repulsive 
nearest-neighbor interactions but does not obey particle- 
hole symmetry (see Fig. 2 in Ref. [22|)- 

It would be interesting to extend the successful treat- 
ment based on the TDFT to higher dimensions. For 
nearest-neighbor interactions there exists a lattice fun- 
damental measure form of the exact zero and one- 
dimensional density functionals [42l | , which enable an ex- 
tension of these functionals to higher dimensions. Al- 
ternatively, the approach used in Sec. IIVI can also be 
generalized to higher dimensions [32| . Based on the re- 
sulting approximate functionals one could make contact 
to previous studies of nearest-neighbor interacting driven 
lattice gases in two dimensions. In these studies struc- 
tural patterns in the NESS were found [sj, as, for ex- 
ample, alternating regions of low and high density for 
attractive interactions V < Ojmanifesting themselves in 
backgammon- [ijl or stripe- [4J| like structures. For re- 
pulsive interactions T^ > 0, the bias can induce transi- 
tions from an ordered to a disordered state [3J, Ea • The 



11 



treatment of these phenomena by TDFT should in partic- 
ular allow one to identify phenomenological parameters 
in former field-theoretical approaches [43, |45| by appro- 
priate coarse-graining. 
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Appendix: Derivation of rate conditions in Eq. (6) 

The master equation 

^^^ = J2 [r(n' ^ n)P(n',i) - r(n ^ n')P(n,t)] 



dt 



(A.l) 



describes the change of the probability P(n, t) of finding 
state n at time t due to transitions from and to other 
states n' with rates r(n' — > n) and r(n —>■ n'), respec- 
tively. For the model in Sec. [IT] with nearest-neighbor 
hopping in the presence of nearest-neighbor interactions, 
we can write 



r(n^ n') 



/l^ 



nih^+iS^, „(i,,+i)r(ni_i, ni+2) , (A. 2) 



where n'^''*+^) is identical to the microstate configuration 
n except that the occupation numbers Ui and ni_|_i are 
interchanged. Inserting this into Eq. (jA.ip . the master 
equation for the stationary state of the TASEP reads 

=^ [n,n,+iF(n(*^*+i)) - n,n,+iF(n)] r(n,_i,n,+2) • 

(A.3) 



Assuming that P{n) ex exp [— 'H(n)] with H from Eq. ([ij, 
P(n(».»+i))/p(n) =exp[-{n^+2-n^^i)V], and Eq. (IX3ll 
becomes 



= > ■ |n,;n,+ie-("'+^-"-i)^ -n,n,+i 



i 



r(ni_i,ni+2) 



0010 ^0100 



r(o,o) 



^oioe-^-iVnooJ r(l,0) 
<oiie''-<ioil r(0,l) 



iV- 



1011 



N, 



101 



r(i,i), 



(A.4) 



where Nq^qq — J2i fii-inini+ihi+2 is the frequency of the 
sequence {0100} of occupation numbers in the microstate 
n, and analogous definitions apply for the remaining N . 
Replacing all n^ by n^ = 1 — n^, the eight numbers N 
can be expressed in terms of the six irreducible num- 
bers iV"x = Y.i^i-T~^ii ^i".i = Z]j"-i-l"»+l' ^ni = 

Y.i'ri'i-inini+i, 7Vi_ii == J^j '^i-i'^J-i-i'^i-l-21 ^11.1 = 
J2i'n'i-inini+2, and N'^m = X]i"»-i'^i'^i-l-i"-i-l-2- This 
yields 

0= [r(i,o)-r(o,i)e^]ivri (a.5) 

+ [r(o,i)-r(i,o)e-^]<i 
+ [r(i, o)e-^ + r(o, i)e^ - r(i, o) - r(o, i)] n^,, 
+ [r(i, G)e-^ + r(o, i)e^ - r(G, 0) - r(i, i)] n^,, 
+ [r(o, 0) + r(i, 1) - r(i, o) - r(o, i)] n;,_, 
+ [r(i, 0) - r(i, o)e-^ + r(o, i) - r(o, i)e^] n^,,, . 

This equation is indeed satisfied for each configuration n 
if the rates fulfill Eqs. (|6ap (vanishing of the first three 
lines and the last line) and (j6bp (vanishing of the fourth 
and fifth line). It is straightforward to extend the analysis 
to ASEPs with detailed balanced backward jump rates 
against the bias direction. Eqs. ((5a|) . (|6b| then specify 
the conditions for the forward rates. 
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